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Introduction 


The  overall  goal  of  this  multi-year  research  project  in  collaboration  with  the  Walter  Reed  Army  Medical  Center 
is  to  develop  the  necessary  technology  to  make  the  proton  facility  that  is  being  constructed  in  Philadelphia  the 
most  advanced  proton  radiotherapy  center.  The  first  technology  is  the  development  of  a  multileaf  collimator 
(MLC)  for  proton  therapy  and  investigates  the  issues  that  must  be  resolved  to  use  an  MLC  in  proton  therapy. 
The  second  technology  under  study  is  the  optimization  of  the  spot-scanning  delivery  technique  including  the 
effects  of  organ  motion.  The  third  technology  is  the  development  of  protocols  to  apply  the  techniques  of  image- 
guided  and  adaptive  radiotherapy  to  proton  therapy,  and  to  develop  a  decision-making  algorithm  to  maximize 
the  efficiency  of  the  facility.  This  report  describes  the  progress  during  the  sixth  year  of  the  expected  seven-year 
process.  Included  in  that  progress  are  the  following  activities  and  achievements:  (1)  Use  of  the  GEANT4  Monte 
Carlo  code,  which  was  developed  in  the  previous  years  of  the  project,  to  test  various  MLC  designs,  culminating 
in  the  delivery  of  the  first  MLC  and  the  status  of  the  on-going  tests  of  that  MLC;  (2)  Use  of  the  same  simulation 
program  to  optimize  the  dose  distribution  from  scanned  beams  accounting  for  inhomogeneities  and  organ 
motion;  and  (3)  Development  of  treatment  protocols  and  an  understanding  of  the  factors  that  are  involved  to 
efficiently  utilize  the  beam. 
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Body 

In  June  2006,  following  years  of  defining  specifications  and  evaluating  proposals,  the  University 
of  Pennsylvania  Health  System  (UPHS)  signed  a  contract  with  Ion  Beam  Applications,  S.A. 
(IBA).  In  addition  to  the  details  associated  with  the  delivery  of  a  proton  therapy  system  the 
contract  included  three  development  agreements  directly  related  to  the  work  supported  by  this 
grant  to  develop  technology  for  proton  therapy.  The  development  agreements  between  UPHS, 
IBA  and  Varian  Medical  Systems,  Inc.  (the  leading  conventional  radiotherapy  vendor)  were:  (1) 
to  develop  a  multileaf  collimator  for  the  IBA  proton  delivery  system,  (2)  to  develop  a  cone-beam 
CT  to  permit  imaging  of  the  patient  in  the  treatment  room,  and  (3)  to  develop  the  pencil-beam 
scanning  algorithm  of  the  Varian  treatment  planning  system. 


Much  of  the  effort  in  the  past  year  has  been  to:  (1)  test  the  prototype  MLC  and  then  install  and 
test  the  final  MLC;  (2)  design  a  system  that  permits  the  treatment  of  shallow  targets  with  scanned 
beams;  and  (3)  write  treatment  protocols  and  submit  them  to  regulatory  bodies  for  approval.  To 
that  end  UPHS  personnel  have  met  with  IBA  and  Varian  engineers  multiple  times  and  have 
teleconferences  or  WebEx-type  remote  meetings  nearly  every  week.  The  MLC,  which  has  the 
highest  priority  because  the  treatment  rooms  cannot  be  commissioned  without  it,  is  the  most 
advanced  of  these  projects.  The  scanning  development  is  close  to  completion;  the  clinical 
implementation  had  been  delayed  by  approximately  one  year  because  we  did  not  feel  that  the 
IBA  system  was  mature  enough  for  routine  clinical  use.  The  cone-beam  CT  development  has 
made  the  least  progress  thus  far,  but  we  continue  to  work  with  several  vendors  and  other  proton 
facilities  to  find  an  affordable  solution.  Cone-beam  CT  technology  has  only  recently  been 
introduced  to  conventional  radiotherapy  and  is  constantly  being  upgraded.  Our  challenge  is  to 
design  a  device  that  will  be  able  to  easily  follow  the  advances  the  system  makes  in  conventional 
therapy. 

This  report  concentrates  on  the  sixth  year  achievements  of  the  multileaf  collimator  development, 
the  fifth  year  of  work  on  the  spot-scanning/motion  project,  and  the  fourth  year  on  the 
development  of  image-guided  and  adaptive  radiotherapy  protocols.  The  Statement  of  Work  in 
the  approved  grant  proposals  included  the  following  items  to  be  investigated.  (Note:  to  minimize 
confusion,  the  years  in  which  we  expected  to  perform  the  work  have  been  replaced  by  the  fiscal 
year  because  there  are  three  separate  starting  dates.)  Because  of  the  delay  in  choosing  the  vendor 
several  of  the  aims  that  were  originally  planned  to  be  completed  by  now,  especially  related  to  the 
design  and  prototyping  of  the  MLC  system,  are  still  ongoing.  The  current  schedule  is  such  that 
the  fifth  treatment  room  (Proton  2  with  pencil  beam  scanning)  will  be  commissioned  in  late 
2011.  The  items  in  the  Statement  of  Work  are  listed  below  with  a  comment  on  the  status  of  any 
item  that  was  to  have  work  performed  by  this  time. 


The  projects  identified  included: 

1 .  Multileaf  collimator  (MLC)  for  use  on  proton  therapy  gantries 

2.  Cone  Beam  CT  (CBCT)  on  the  proton  gantry  for  localization  of  target  volumes 

3.  Proton  Radiography  to  determine  dose  and  stopping  power  of  various  tissues 

4.  Positron  Emission  Tomography  (PET)  imaging  on  the  gantry  to  evaluate  dose  deposition 
within  tissues  irradiated 
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5.  Pencil  Beam  Scanning  (PBS)  proton  beam  using  adaptive  radiotherapy  techniques  based 
on  implementation  of  MLC,  Cone  Beam  CT,  and  PET  imaging 


Phase  One-  All  the  action  items  in  the  SOW  for  phase  one  were  completed  by  us  with  minimal 
adjustments  in  some  of  the  items  that  were  required  by  the  development  of  the  proton  system. 

Phase  Two-  The  SOW  for  phase  two  included  mostly  the  development  and  implementation  of 
Pencil  Beam  Scanning,  and  development  of  joint  military  UPENN  telemedicine  system.  The 
delivery  and  installation  of  clinical  PBS  system  by  our  contractor-  IBA-  has  been  delayed  by 
IBA.  The  first  PBS  system  is  currently  being  installed  and  tested  at  our  facility  by  IBA.  All  the 
action  items  in  the  SOW  that  are  related  to  our  preparation  for  PBS  and  could  be  done  at  the 
current  stage  of  installation  have  been  done.  The  telemedicine  system  with  WRMC  is  working. 
We  expect  to  have  our  first  clinical  PBS  system  before  the  end  of  the  year  and  then  will  be  able 
to  complete  the  action  items  that  are  in  the  SOW  of  this  phase. 

Phase  Three-  The  SOW  included  studying  of  imaging  modalities  for  proton,  development 
of  CBCT,  studying  and  development  of  adaptive  proton  therapy  which  is  based  on  repeated 
imaging  during  proton  therapy  and  development  of  smart  scheduling  system. 

The  studies  of  imaging  modalities  per  the  SOW  were  done. 

We  currently  do  not  have  CBCT  despite  intensive  work  and  research  to  develop  one.  We  are 
actively  reviewing  new  idea  to  develop  CBCT  with  IBA  and  have  successfully  completed 
preliminary  testing.  In  this  proposal  the  CBCT  will  be  installed  on  the  gantry.  We  plan  to 
continue  and  develop  this  idea  with  IBA. 

To  satisfy  the  SOW  regarding  adaptive  radiotherapy  for  proton  we  have  developed  IRB 
approved  protocol  which  is  currently  being  reviewed  by  the  DOD.  Once  approved  we  will  recruit 
the  additional  personal  needed  to  execute  it. 

Current  and  near  future  plan  to  complete  the  SOW: 

Complete  the  items  that  are  related  to  PBS  after  the  expected  installation  of  PBS  during  the 
second  quarter  of  2012,  as  was  detailed  above. 

We  continue  to  review  the  new  proposal  for  cone  beam  CT  and  we  are  expecting  to  proceed  with 
the  development  of  cone  beam  CT  during  the  first  half  of  FY2012. 

Once  final  DOD  approval  is  received  we  will  continue  the  development  of  adaptive  radiotherapy 
methodology  which  will  continue  until  September  of  2012. 

We  continue  to  develop  the  smart  scheduling  system  as  the  department  grows  and  plan  to 
complete  it  by  September  of  2012 


Progress 

The  work  over  the  last  year  can  be  classified  into  three  areas  relating  to  progresses  in:  (A)  MLC 
development,  (B)  spot-scanning  development,  and  (C)  protocol  development. 


A.  MLC  progress 
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I.  Small  field  measurements 


In  our  previous  report,  we  compared  simulation  results  with  small  field  measurements  for 
a  pristine  peak  with  a  range  of  24.67  cm  using  circular  MLC  shapes  with  several  different 
diameters.  The  measurements  were  taken  with  a  PPC05  ion  chamber  and  the  depth  dose 
curves  from  simulation  were  found  by  integrating  through  a  lxl  cm2  square  in  the  center 
of  the  water  tank.  We  also  made  comparisons  with  measurements  and  Eclipse 
calculations  of  SOBPs  with  a  range  of  17.5  cm  and  a  10  cm  modulation  for  a  number  of 
circular  MLC  shapes  differing  in  diameter.  The  measured  SOBP  for  a  1.8  cm  diameter 
circular  field  size  differed  significantly  from  the  Eclipse-calculated  SOBP  for  the  same 
field  size.  Measurements  showed  a  larger  decrease  in  dose  for  the  small  field  size  in 
comparison  with  the  SOBPs  for  larger  field  sizes,  although  both  Eclipse  and 
measurement  showed  that  the  SOBP  for  the  1 .8  cm  diameter  circular  field  had  a 
significant  downward  slope  in  comparison  with  the  SOBPs  for  larger  field  sizes.  We 
found  that  as  we  decreased  the  area  of  the  center  of  the  water  phantom  through  which  the 
simulation  is  integrated  below  lxl  cm2,  the  simulated  SOBP  deviated  from  the  chamber- 
measured  SOBP  and  approached  the  Eclipse-calculated  SOBP. 

In  order  to  better  understand  the  tilt  in  the  SOBP  for  small  field  sizes,  we  looked  at  x-  and 
y-profiles  and  XY-dose  distributions  at  multiple  depths  for  the  1.8  cm  diameter  circular 
field  and  the  5  cm  diameter  field  (for  comparison).  The  downward  tilt  implies  that  the 
long  range  protons  are  not  depositing  dose  along  the  central  axis,  so  we  want  to  look  at 
where  they  are  traveling.  Figure  1.1  contains  one  of  the  SOBP  comparisons  described 
above  from  the  last  report  containing  measured  and  Eclipse-calculated  SOBPs  for  various 
field  sizes  with  the  simulated  SOBPs  found  using  a  6x6mm2  integration  area  through  the 
center  of  the  water  phantom.  Figure  1.2  contains  the  XZ  and  YZ,  2D  projections  of  the 
dose  through  the  phantom  which  have  been  integrated  through  the  y-  and  x-direction, 
respectively.  The  XY-dose  distributions  and  profiles  at  different  depths  for  the  beams 
from  Figure  1.1  are  shown  in  Figure  1.3.  The  XY-dose  distributions  have  a  6x6mm2  square 
in  black  in  the  center  to  give  an  idea  of  the  area  being  integrated  through  the  z-direction 
to  produce  the  simulated  depth  dose  curves.  The  ID  profiles  in  the  y-direction  have  the 
Eclipse-calculated  y-profile  at  isocenter  (d  =  12.5  cm)  plotted  for  comparison  of  shape. 

In  addition  to  these  small  field  measurements,  we  simulated  a  9.55  cm  range  pristine  peak 
with  0.86  cm  side  square  field  size  for  which  measurements  had  been  taken  with  a  PPC05 
ion  chamber  and  an  SFD  stereotactic  diode.  The  chamber  and  diode  measurements 
differed  by  approximately  20%  in  entrance  dose,  with  the  chamber  dose  being  lower  than 
the  diode.  The  measured  peaks  were  plotted  with  the  simulated  pristine  peak  obtained 
using  12x1 2mm  ,  1 0x1 0mm  ,  8x8mm  ,  6x6mm  ,  and  4x4mm  integration  areas  through 
the  center  of  the  water  phantom.  We  found  that  the  simulated  peak  matched  with  the 
chamber  measurement  when  using  an  integration  area  of  12x12  mm  through  the  center 
of  the  water  phantom  (which  is  a  greater  area  than  the  field  size),  and  that  after  shrinking 
the  integration  area  to  4x4mm2,  the  simulation  entrance  dose  matched  with  that  of  the 
pristine  peak  measured  using  a  diode.  Although  using  an  integration  area  smaller  than 
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4x4mm2  was  attempted,  the  data  was  too  noisy  and  would  require  smaller  voxels  to  be 
used  in  the  phantom  to  get  useful  statistics.  The  voxels  that  were  used  in  this  simulation 
were  0.5x0. 5x1.0  mm3  in  a  20x20x40  cm3  water  phantom.  The  results  can  be  seen  in 
Figure  1.4  which  shows  the  peaks  measured  with  the  chamber  and  diode  as  well  as  the 
simulated  peak  as  a  function  of  integration  area  through  the  center  of  the  water  tank.  This 
shows  that  the  chamber  may  not  produce  accurate  measurements  for  small  field  sizes. 

B6  range  =  17.5cm,  modulation  =  10cm,  0.6cm  integration 


Figure  AI.l:  A  comparison  between  SOBPs  for  multiple  circular  field  sizes  measured  in  a  water  tank 
with  a  PPC05  chamber  (dotted  lines),  calculated  in  Eclipse  (dashed  lines),  and  simulated  in  a  water 
phantom  (solid  lines).  The  SOBPs  have  a  range  of  17.5  cm  and  a  modulation  of  10  cm  and  were 
measured  for  circular  field  sizes  with  10  (purple),  3  (blue),  and  1.8  (blue)  cm  diameters  and  simulated 
using  circular  field  sizes  of  5,  3,  1.8,  and  1  cm  diameters.  There  is  no  noticeable  difference  for  field 
sizes  with  diameters  greater  than  5  cm,  so  the  simulation  with  5  cm  diameter  field  size  corresponds 
with  the  measured  SOBP  for  a  10  cm  diameter  field  size.  The  SOBPs  from  the  simulations  have  been 
obtained  by  integrating  the  dose  over  a  6x6  mm2  square  through  the  center  of  the  water  phantom. 

I. 
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Figure  AI.2:  2-Dimensional  dose  distributions  for  an  SOBP  with  a  range  of  17.5  cm,  modulation  of  10 
cm,  and  circular  field  sizes  with  a  1.8  cm  (top)  and  5  cm  (bottom)  diameter.  The  XZ-projections  (left) 
have  been  integrated  through  the  y-direction  and  the  YZ-projections  (right)  have  been  integrated 
through  the  x-direction.  They  have  been  normalized  to  the  dose  in  the  center  of  the  SOBP  and  the  dose 
has  been  integrated  through  the  entire  phantom  in  the  respective  directions.  They  give  a  qualitative  idea 
of  where  the  long  range  protons  are  traveling  for  the  small  field  in  comparison  with  the  larger  field. 
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Figure  A1.3:  XY-dose  distributions  and  profiles  in  the  x-  and  y-  direction  for  circular  fields  with 
diameters  of  1.8  cm  (left)  and  5  cm  (right)  at  depths  of  0,  5,  10,  12.5,  15,  and  17  cm.  The  beam  is 
double  scattering  option  B6  with  a  range  of  17.5  cm  and  a  10  cm  modulation.  There  is  a  square  in  the 
center  of  the  XY-dose  distributions  indicates  the  area  of  integration  used  to  create  the  simulated  SOBPs 
in  Figure  1.  The  y-profiles  have  the  Eclipse-calculated  profile  at  isocenter  (depth  =  12.5  cm)  plotted  for 
reference  in  purple. 
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range  =  9,55  cm,  12  mm  integration 


range  =  9,55  cm,  10  mm  integration 


range  =  9.55  cm,  8  mm  integration 


range  =  9.55  cm,  6  mm  integration 
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ninge  =  9.55  cm,  4  mm  integration 


Figure  AI.4:  Pristine  peaks  using  double  scattering  option  B6  with  a  range  of  9.55  cm  and  a 
0.86x0. 86cm2  square  field  size  comparing  simulation  curves  (purple)  with  measurements  taken  with 
both  a  PPC05  ion  chamber  (green)  and  an  SFD  stereotactic  diode  (red).  The  simulated  pristine  peaks 
have  been  plotted  using  different  areas  of  integration  through  the  center  of  the  water  phantom.  The 
areas  that  have  been  used  are  12x12,  10x10,  8x8,  6x6,  and  4x4  mm2.  With  an  integration  area 
~12xl  2mm2  the  simulation  is  comparable  to  the  chamber  measurement,  but  as  the  integration  area 
decreases,  the  simulation  gets  closer  to  the  diode  measurement. 

II.  Simulating  the  compensator 

At  the  time  of  our  last  report,  we  were  having  trouble  creating  patient-specific 
compensators  in  the  simulation  because  of  the  complicated  array  of  volume  boundaries 
created  when  using  a  large  number  of  subtraction  volumes.  Geant  has  trouble  handling 
the  tracking  of  particles  through  overlapping  Boolean  volumes,  and  because 
compensators  have  thousands  of  drill  positions,  our  simulations  using  subtraction 
volumes  to  create  the  compensator  were  unable  to  finish  after  several  days.  In  order  to 
solve  this  problem,  we  decided  to  try  to  approximate  the  compensator  such  that  we  did 
not  use  overlapping  or  subtraction  volumes.  We  did  this  by  using  the  remaining  thickness 
of  the  compensator  at  each  dwell  position  to  create  a  rectangular  prism  that  did  not 
overlap  with  neighboring  dwell  positions.  A  diagram  of  the  pattern  used  is  shown  in 
Figure  II.  1 .  Because  clinical  compensators  are  constructed  with  several  thousand  drill 
positions,  this  approximation  should  not  cause  much  difference  in  dose  distribution. 

We  tested  the  orientation  of  the  compensator  by  creating  a  test  compensator  with  four 
different  thicknesses  in  each  comer.  The  test  compensator  was  2mm  thick  in  the  upper 
left  comer,  20mm  in  the  upper  right  comer,  40mm  in  the  lower  right  comer,  and  60mm  in 
the  lower  left  comer.  The  orientation  of  the  compensator  with  respect  to  the  MLC  was 
verified  by  visualizing  and  measuring  the  thicknesses  of  the  compensator  using  HepRep. 
Figure  II.2  shows  the  XY-dose  distributions  at  different  depths  for  the  test  compensator, 
again  verifying  the  orientation  of  the  compensator  created  in  simulation.  The  file  used  to 
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create  this  compensator  had  over  17,500  drill  positions,  meaning  there  were  over  17,500 
volumes  created  in  the  process  of  making  the  simulated  compensator.  In  this  case,  where 
we  are  no  longer  using  overlapping  subtraction  volumes  to  make  the  compensator,  the 
simulation  did  not  take  much  longer  to  run  than  if  the  same  simulation  had  been  run  with 
the  compensator  absent.  The  simulation  can  now  read  the  drill  positions  and  remaining 
Plexiglas  thicknesses  of  the  compensator  from  a  DICOM  file,  and  create  the  field- 
specific  compensators  in  our  model. 

When  running  patient  simulations,  we  found  that  if  the  air-gap  between  the  compensator 
and  the  patient  was  small,  as  for  brain  patients,  when  the  CT  images  were  reconstructed 
in  simulation  the  field  of  view  extended  such  that  the  compensator  was  being  overwritten 
by  air  material  due  to  the  large  area  of  air  in  the  CT  set  surrounding  the  patient.  For  this 
reason  a  crop  function  was  developed  so  that  the  CT  set  can  be  cropped  around  the 
patient  and  the  compensator  will  no  longer  be  overwritten. 
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Figure  AII.l:  To  create  the  compensator,  the  simulation  reads  in  the  row  and  column  spacing  and  row 
offset  of  the  drill  positions  from  a  DICOM  file  and  creates  a  rectangular  prism  with  the  remaining 
thickness  of  Plexiglas  at  each  position.  Previously  we  had  been  creating  the  compensator  by  making 
cylindrical  subtraction  volumes  with  the  milling  tool  diameter,  but  Geant  was  having  problems 
handling  the  large  number  of  overlapping  Boolean  volumes. 

Test  Pattern  Depth  =  20  cm 
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Figure  AII.2:  To  test  that  the  compensator  is  being  created  in  the  correct  orientation,  a  test 
compensator  was  created  following  the  pattern  shown  in  the  upper  left  with  different  thicknesses  of 
Plexiglas  at  each  comer  using  17,555  drill  positions.  A  beam  with  a  range  of  25  cm  was  sent  through 
the  compensator  and  the  XY-dose  distribution  at  depths  of  20,  22,  and  24  cm  are  shown  here.  At  a 
depth  of  20  cm  in  the  comer  where  the  compensator  is  60mm  there  is  almost  no  dose,  and  in  the  comer 
where  the  compensator  is  40mm  thick,  the  peak  dose  is  being  deposited.  At  a  depth  of  22cm,  the  comer 
corresponding  to  where  the  compensator  is  40mm  thick  has  greatly  reduced  dose,  and  in  the  comer 
where  is  compensator  is  20mm  thick,  the  peak  dose  is  being  deposited.  This  shows  that  the 
compensator  is  being  created  in  simulation  with  the  correct  orientation. 

III.  Optimization  of  the  Simulation  Model 

The  last  few  quarterly  reports  have  been  focused  on  the  verification  of  our  simulation 
model  with  measured  beam  data  so  that  simulations  can  be  used  in  the  future  for  clinical 
investigations  with  confidence  in  the  accuracy  of  the  results.  At  this  point  we  have 
double-scattering  beam  data  from  three  different  treatment  rooms  as  well  as  uniform 
scanning  data  from  one  room  that  we  are  utilizing  to  investigate  differences  between 
simulation  results  and  measured  data  as  well  as  to  tune  the  simulation  beam  parameters 
such  that  the  results  produced  more  closely  resemble  measured  data. 

In  our  previous  report,  we  discussed  changing  from  Geant  version  4.9.1  to  version  4.9.3. 

Due  to  continued  discrepancies  between  the  range  of  the  simulated  beam  and  measured 
data,  we  have  updated  to  the  newest  version  of  Geant  (4.9.4).  We  had  found  that 
simulated  results  from  version  4.9.3  closely  matched  the  range  of  uniform  scanning 
pristine  peaks  that  were  measured  without  anything  in  the  beam  line,  indicating  that  it 
accurately  modeled  the  beam  range  in  water.  However,  once  the  double-scattering  beam 
line  components  were  added  in  the  simulated  beam  line,  the  resulting  simulated  beam 
range  was  consistently  deeper  than  the  measured  range.  In  changing  the  Geant  version 
from  4.9.3  to  4.9.4  we  wanted  to  first,  verify  that  range  of  the  simulated  beam  in  water 
without  any  beam  line  components  did  not  change,  then  check  if  the  range  of  double- 
scattered  pristine  peaks  were  any  closer  to  measured  ranges. 
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As  described  in  our  previous  quarterly  report,  pristine  peaks  measured  during  uniform 
scanning  commissioning  were  compared  with  simulation  results  employing  nominal 
spectral  data  provided  by  the  manufacturer  and  the  same  beamline  settings.  The  measured 
uniform  scanning  peaks  with  nothing  in  the  beam  line  were  used  to  tune  the  energy 
spectrum  of  the  beam  in  the  Monte  Carlo  such  that  residual  differences  between 
simulation  and  measurement  were  minimized.  Because  the  entrance  dose  in  simulation 
was  found  initially  to  be  consistently  lower  than  in  measurement,  a  distribution  of  low 
energy  protons  as  suggested  by  Bortfeld  (MedPhys  1997)  was  added  to  the  beam  to 
represent  protons  that  may  have  scattered  from  the  edges  of  the  beam  pipe,  but  would  not 
otherwise  be  accounted  for  in  the  simulation  model. 

Bortfeld’s  theory  was  that  the  beam  can  be  split  into  a  narrow  Gaussian  distribution  about 
the  mean  kinetic  energy  of  the  beam  and  a  tail  of  low  energy  protons  due  to  scatter  along 
the  beam  pipe  which  is  a  small  percentage  of  the  protons  in  the  Gaussian.  Because  the 
fluence  in  the  low  energy  tail  is  only  a  small  fraction  of  the  total  beam  fluence  and  it  is 
expected  that  the  fluence  of  protons  with  zero  energy  will  be  zero  and  there  will  be  a 
positive  ramp  of  fluence  up  to  the  mean  kinetic  energy,  the  low  energy  tail  can  be 
approximated  using  a  linear  ramp  such  that  <DE  (E)  a  E  for  0<E<Eo.  To  better  fit  the 
widths  of  the  pristine  peaks,  the  sigma  of  the  Gaussian  of  the  beam’s  mean  kinetic  energy 
was  also  adjusted.  We  found  that  by  finding  the  best  fitting  kinetic  energy  spread  and 
then  varying  the  fraction  the  low  energy  proton  tail  to  be  added,  we  could  get  the 
simulated  pristine  peaks  to  fit  perfectly  with  the  measured  data. 

The  spread  of  the  Gaussian  about  the  mean  kinetic  energy  had  previously  been  set  as  a 
function  of  percent  momentum  spread  which  was  calculated  as  follows: 

Percent  momentum  spread  =  0.00801  -  3.08E-5*Mean  kinetic  energy 

where  the  two  parameters  for  this  linear  fit  come  from  beam  data  originally  provided  by 
IBA.  Several  simulations  were  run  corresponding  to  each  of  the  uniform  scanning 
measurements  with  nothing  in  the  beam  line  varying  the  kinetic  energy  spread  to  find  the 
optimal  fit.  The  dose  distribution  of  the  low  energy  tail  was  then  added  to  the  results  of 
the  beam  with  the  optimal  kinetic  energy  spread  to  optimize  the  fit  of  simulation  results 
to  measurement. 

To  verify  that  the  range  in  water  did  not  change,  we  repeated  simulations  that 
corresponded  to  9  pristine  peaks  taken  during  uniform  scanning  commissioning  with 
nothing  in  the  beamline.  The  results  from  simulations  with  version  4.9.4  were  plotted 
along  with  measured  pristine  peaks  and  the  results  from  simulations  using  Geant  version 
4.9.3.  The  results  are  plotted  in  Figure  III.  1  and  they  show  that  there  was  no  significant 
change  in  the  range  of  protons  in  water  between  the  two  versions. 

Next  we  wanted  to  check  how  the  double-scattering  SOBPs  using  Geant  4.9.4  looked 
with  respect  to  the  measurements  from  the  three  rooms  that  have  been  commissioned,  so 
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simulations  were  run  for  a  low,  medium,  and  high  energy  range  SOBP.  Since  the  beam 
current  modulation  files  that  we  are  using  were  optimized  to  create  flat  SOBPs  when  we 
were  using  Geant  4.9.1,  we  did  not  expect  the  simulated  SOBPs  using  Geant  4.9.4  to  be 
flat.  The  results  are  plotted  in  Figure  III.2  and  the  same  results  using  Geant  4.9.3  are 
plotted  for  comparison  to  show  that  the  newer  version  made  an  improvement  concerning 
the  range  of  double-scattering  SOBPs. 

All  of  the  double-scattering  pristine  peaks  that  were  measured  during  commissioning  of 
rooms  1,  3,  and  5  were  also  simulated  using  Geant  4.9.4.  The  simulations  were  run 
incorporating  the  optimized  kinetic  energy  sigma  and  fraction  of  the  low  energy  proton 
distribution  that  were  described  in  the  previous  report.  These  simulations  were  run  to 
check  the  range  of  the  simulated  results  and  to  see  how  the  fit  for  the  kinetic  energy 
sigma  and  low  energy  proton  distribution,  which  were  optimized  to  fit  uniform  scanning 
pristine  peaks,  affected  the  double-scattering  simulations  results.  The  simulations  results 
are  plotted  in  Figure  III.3  with  the  measured  pristine  peaks  for  each  of  the  three  treatment 
rooms  for  comparison.  Table  III.  1  lists  the  nominal  distal  90%  range,  the  average  range 
from  commissioning  measurements  in  three  treatment  rooms,  and  the  ranges  from 
simulation.  It  also  lists  the  difference  between  the  range  from  the  simulation  results 
(versions  4.9.1  and  4.9.4)  and  measured  range,  and  the  differences  are  plotted  as  a 
function  of  nominal  range  in  Figure  III.4. 

For  each  of  the  Geant  versions  that  we  have  tested  for  use  in  our  simulation  (  4.9.1,  4.9.3, 
and  4.9.4  )  the  water-equivalent  thickness  (WET)  of  the  second  scatterers  and  modulator 
wheel  have  been  found  by  running  simulations  for  each  of  the  double-scattering  options 
first  with  nothing  in  the  beam  line,  and  then  with  each  of  the  components  being  placed  in 
the  beam  line.  The  WET  of  the  components  is  of  interest  because,  as  mentioned 
previously,  the  simulated  range  is  equal  to  the  measured  range  in  water  when  there  is 
nothing  in  the  beam  line,  but  when  looking  at  double-scattering  measurements,  the 
simulated  range  is  too  long.  This  means  that  the  range  of  the  beam  is  not  being  reduced 
by  the  same  amount  in  simulation  as  in  reality  while  traveling  through  the  beam  line 
components.  Table  III.2  contains  the  WET  of  the  second  scatterers  and  modulator  wheels 
given  by  IB  A  and  the  WETs  found  in  simulation  for  each  double-scattering  energy 
options  for  each  of  the  Geant  versions  that  we  have  tested. 


19 


US2  range  =  14  cm 


US2  range  =  16  cm 


Ccniir  9.3 
-  M  cm  mi  i  v  mvii  ■ 


Z  /dll 


US2  range  =  18  cm 


IJS3  range  =  20-5  cm 


20 


US3  range  =  23  cm 


US3  range  =  25.99  cm 


US4  range  =  28.5  cm 


US4  range  —  3ft.5  cm 


21 


US4  range  =  32.14  cm 


Figure  AIII.l:  Uniform  scanning  pristine  peaks  simulated  with  Geant  4.9.3  (purple),  Geant  4.9.4 
(blue)  plotted  with  measured  pristine  peaks  (black)  to  verify  that  the  range  in  water  between  the  two 

versions  does  not  change.  Measurements  and  simulation  were  done  with  a  lOxlOcm^  square  field  size 

with  a  40x40x40cm3  water  tank.  For  these  particular  ranges  in  uniform  scanning  there  are  no 
scatterers  in  the  beam  line  and  the  beam  does  not  travel  through  the  modulator  wheel.  The  simulations 
were  run  using  the  newly  optimized  kinetic  energy  spread  mentioned  in  the  previous  report. 
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Figure  AIII.2:  Double-scattering  depth  dose  curves  for  SOBPs  that  have  been  measured  during 
commissioning  for  rooms  1,  3,  and  5  plotted  with  results  from  simulations  using  Geant  4.9.3 
(left)  and  Geant  4.9.4  (right).  The  simulation  results  are  plotted  in  black. 
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Figure  AIII.3:  Double-scattering  pristine  peaks  that  have  been  simulated  using  Geant  4.9.4  with  the 
newly  optimized  sigma  of  the  kinetic  energy  Gaussian  about  the  mean  kinetic  energy  and  the  added 
distribution  of  low  energy  protons.  Simulation  results  are  plotted  in  red  while  the  commissioning 
measurements  from  the  three  treatment  rooms  are  plotted  in  black. 


Table  AIIL1:  Summary  of  ranges  from  double-scattering  pristine  peaks 


Nominal  Range 

Avg  Measured(Pl,4,5) 

G9.1 

G9.4 

G9.1-MeasuredAvg 

G9.4-MeasuredAvg 

B1 

4.6 

4.647 

5.1 

4.9 

0.453 

0.253 

4.99 

5.040 

5.5 

5.3 

0.460 

0.260 

5.41 

5.450 

5.9 

5.7 

0.450 

0.250 

5.86 

5.910 

6.3 

6.2 

0.390 

0.290 

B2 

5.87 

5.933 

6.3 

6.2 

0.367 

0.267 

6.37 

6.423 

6.9 

6.7 

0.477 

0.277 

6.91 

6.957 

7.4 

7.2 

0.443 

0.243 

7.49 

7.553 

8 

7.8 

0.447 

0.247 

B3 

7.5 

7.587 

7.9 

7.8 

0.313 

0.213 

8.12 

8.183 

8.5 

8.4 

0.317 

0.217 

8.81 

8.887 

9.2 

9.1 

0.313 

0.213 

9.54 

9.620 

9.9 

9.8 

0.280 

0.180 

B4 

9.55 

9.540 

10 

9.9 

0.460 

0.360 

10.32 

10.443 

10.8 

10.7 

0.357 

0.257 

11.09 

11.217 

11.6 

11.4 

0.383 

0.183 

11.85 

11.950 

12.3 

12.2 

0.350 

0.250 

B5 

11.86 

12.030 

12.2 

12.3 

0.170 

0.270 

26 


13.09 

13.243 

13.4 

13.5 

0.157 

0.257 

14.32 

14.477 

14.6 

14.7 

0.123 

0.223 

15.53 

15.677 

15.9 

15.9 

0.223 

0.223 

B6 

15.54 

15.640 

15.8 

15.9 

0.160 

0.260 

16.86 

16.977 

17.1 

17.2 

0.123 

0.223 

18.28 

18.383 

18.5 

18.6 

0.117 

0.217 

19.83 

19.960 

20.1 

20.1 

0.140 

0.140 

B7 

19.84 

20.017 

20.1 

20.2 

0.083 

0.183 

21.18 

21.350 

21.5 

21.5 

0.150 

0.150 

22.61 

22.690 

23.1 

22.9 

0.410 

0.210 

23.91 

24.020 

24.1 

24.2 

0.080 

0.180 

B8 

22.8 

22.873 

22.9 

22.9 

0.027 

0.027 

24.67 

24.730 

24.7 

24.7 

-0.030 

-0.030 

26.54 

26.590 

26.6 

26.6 

0.010 

0.010 

28.26 

28.237 

28.3 

28.3 

0.063 

0.063 

V. 


Difference  Between  Simulation  and  Measurement  vs  Nominal 

Range 


0.600 

0.500 

0.400 

0.300 

0.200 


2  o.ioo 

3 
£ 

0.000 

-0.100 


♦♦♦ 


♦  ♦ 

♦  ♦  ♦ 


10 


15 


20 


?5 


30 


Nominal  Range  (cm) 


♦  G9.1-MeasuredAvg 
■  G9.4 — MeasuredAvg 


Table  AIII.2:  List  of  simulated  and  IBA-provided  WET  of  the  second  scatterers 
and  modulator  wheels  for  each  of  the  double-scattering  options  using  each  Geant 
version  we  have  tested. 


Geant  4.9.1 


MC  WET 

ss 

IBA  WET 

SS 

Difference(cm) 

MC  WET 
ModWheel 

IBA  WET 
ModWheel 

Difference 

(cm) 

B1 

6.1 

6.2000 

-0.1 

4.7 

4.83 

-0.13 

B2 

6.2 

6.2500 

-0.05 

4.7 

4.83 

-0.13 
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B3 

6.1 

6.2500 

-0.15 

1.7 

1.78 

-0.08 

B4 

6.2 

6.2800 

-0.08 

2.2 

2.23 

-0.03 

B5 

3.5 

3.6000 

-0.1 

1.7 

1.78 

-0.08 

B6 

3.5 

3.6100 

-0.11 

2.2 

2.23 

-0.03 

B7 

3.6 

3.6400 

-0.04 

2.9 

2.83 

0.07 

B8 

0.8 

0.8400 

-0.04 

1.9 

1.89 

0.01 

Geant  4.9.3 


MC  WET 

ss 

IBA  WET 
SS 

Difference( 

cm) 

MC  WET 
ModWheel 

IBA  WET 
ModWheel 

Difference( 

cm) 

B 

1 

6.1000 

6.2000 

-0.1 

4.6 

4.83 

-0.23 

B 

2 

6.1000 

6.2500 

-0.15 

4.6 

4.83 

-0.23 

B 

3 

6.1000 

6.2500 

-0.15 

1.6 

1.78 

-0.18 

B 

4 

6.1000 

6.2800 

-0.18 

2 

2.23 

-0.23 

B 

5 

3.3000 

3.6000 

-0.3 

1.6 

1.78 

-0.18 

B 

6 

3.3000 

3.6100 

-0.31 

2 

2.23 

-0.23 

B 

7 

3.4000 

3.6400 

-0.24 

2.7 

2.83 

-0.13 

B 

8 

0.8000 

0.8400 

-0.04 

1.7 
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B.  Spot  Scanning  Progress 

I.  Eclipse  Plan  Verification  with  Monte  Carlo 

A  major  accomplishment  of  the  Phase  II  research  this  past  year  is  the  development  of 
software  that,  given  a  patient  CT  dataset  and  a  treatment  plan  file  from  the  Eclipse 
treatment  planning  system,  performs  the  calculation  of  the  3D  dose  distribution  in  patient 
using  the  Geant4  Monte  Carlo  engine.  This  is  called  the  forward  dose  calculation:  given 
the  beam  delivery  specification  (from  the  treatment  plan  file),  the  dose  to  the  patient  is 
determined.  (The  inverse  calc  is  to  determine  the  delivery  parameters  from  the  desired 
dose  distribution  in  patient.)  For  the  case  of  pencil-beam  scanning,  the  Eclipse  treatment 
planning  system  was  commissioned  with  Monte  Carlo  generated  beam  profiles  and  depth 
dose  curves.  The  7.5  cm  acrylic  range  shifter  used  to  treat  shallow  depths  was  also 
simulated  and  commissioned  in  Eclipse.  In  order  to  calibrate  the  number  of  simulated 
protons  per  MU,  a  reference  condition  was  calculated  with  Eclipse:  a  10  cm  x  10  cm  x  10 
cm  target  within  a  water  phantom  was  brought  to  1  Gy  and  the  number  of  protons  per 
MU  was  adjusted  to  make  the  total  number  of  MU  for  this  reference  condition  treatment 
plan  equal  to  1000.  This  is  different  than  the  reference  condition  previously  quoted,  100 
MU  for  1  Gy/L,  that  is  in  use  at  MD  Anderson.  The  particular  value  is  not  important,  as 
long  as  it  is  known  and  stable,  and  that  the  absolute  dose  calibration  factor  (Gy/MU)  is 
measured  accurately.  The  value  of  1000  was  recommended  by  IB  A,  and  the  specs  quoted 
by  IB  A  like  the  minimum  number  of  MU  per  spot  and  the  resolution  of  the  monitor 
chambers  assume  this  definition  for  MU.  The  depth  range  for  the  reference  condition  is 
11  cm  to  21  cm  depth.  For  simulations,  the  number  of  protons  per  MU  should  be  a 
function  of  MU,  because  the  monitor  chamber  response  is  determined  by  the  stopping 
power  of  protons  in  air,  which  is  a  function  of  energy,  but  also  the  mean  path  length  of 
particles  through  the  monitor  chamber  will  be  a  function  of  energy:  the  lower  energy 
beams  have  a  larger  beam  emittance  and  will  tend  to  have  a  slightly  large  path  length. 

The  second  effect  has  not  yet  been  included  in  simulation.  The  first  effect  has  been 
included  by  incorporating  the  stopping  power  in  air  into  the  calculation  of  the  number  of 
protons  per  MU.  A  function  of  the  form  A  +  B/E  *  In  E  is  motiviated  by  taking  the  non- 
relativistic  limit  of  the  intermediate-energy  term  of  the  Bethe-Bloch  equation.  The 
parameters  A  and  B  were  obtained  by  PLSR  fit  to  data  from  NIST  PStar.  It  is  interesting 
that  the  stopping  power  in  air  decreases  by  about  a  factor  of  2  over  the  range  of  100  MeV 
to  230,  and  that  the  height  of  peak  of  the  proton  depth  dose  curve  is  reduced  by  nearly  the 
same  factor  of  2  due  to  range  straggling  in  such  a  way  that  the  two  effects  balance  (at 
least  over  this  energy  range),  and  a  plot  of  Bragg  peaks  of  different  energies  have  roughly 
the  same  Gy/MU  in  the  peak  (see  Figure  1-1). 
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We  can  perform  the  forward  dose  calculation  using  our  in-house  PBS  simulation  code 
based  on  the  Geant4  Monte  Carlo  engine.  We  will  use  these  simulations  to  validate  the 
dose  calc  algorithms  in  Eclipse  and  to  tune  parameters  in  Eclipse  that  have  an  effect  on 
fits  of  the  commissioning  data  and  parameters  affecting  the  optimization  engine.  Derek 
has  commissioned  Eclipse  with  simulated  data  using  the  universal  and  dedicated  nozzle 
codes  so  that  planning  studies  can  proceed  with  Eclipse.  Projects  are  underway  at  Penn  to 
investigate  the  quality  of  PBS  plans,  in  terms  of  how  to  optimize  PBS  plans;  how  plan 
quality  can  be  improved  with  fine  beam  control,  range  shifters  and/or  bolus,  possibly 
lateral  spreaders;  and  how  PBS  plans  will  compare  to  other  modalities  for  specific 
disease  sites  (currently,  prostate,  head  and  neck,  breast,  rectal,  and  pancreatic  cancers). 
An  example  of  the  PBS  forward  calculation  using  Geant  is  shown  in  Figure  1-2.  The 
dose-volume-histograms  (DVH)  for  the  Monte  Carlo  plan  is  much  worse  than  the 
Eclipse-calculated  plan.  This  is  due  to  at  least  two  reasons.  The  Monte  Carlo  calc  needs 
to  be  run  with  more  primaries  for  higher  accuracy.  Dose  is  accurate  to  only  2-5%  in  the 
PTV  (see  Figure  1-3).  The  random  statistical  noise  in  the  Monte  Carlo  dose  distribution 
blurs  the  DVH  line  for  target  coverage.  Still,  the  median  dose  to  the  PTV  is  shifted  in 
Monte  Carlo.  This  is  believed  to  be  the  effect  of  the  dose  halo.  Eclipse  is  modeling  the 
proton  beam  with  a  single  Gaussian,  but  there  is  a  Rutherford  tail  to  the  profile  due  to 
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nuclear  interactions  (E  Pedroni  et  al  2005  Phys.  Med.  Biol.  50  541  doi:  10.1 088/003 1  - 
9155/50/3/011).  Figure  1-4  shows  the  distribution  of  the  difference-of-dose  between 
Monte  Carlo  and  Eclipse.  The  two  calculations  mostly  agree  in  the  center  of  the  PTV. 
Around  the  margin  of  the  PTV,  Monte  Carlo  is  about  5%  lower  than  Eclipse,  and  a  little 
further  outside  the  PTV  Monte  Carlo  is  5%  higher.  This  is  consistent  with  the  idea  that 
some  of  the  Gaussian  peak  has  been  scattered  out  to  a  longer  tail  that  Eclipse  is  not 
modeling.  Within  Eclipse,  it  is  possible  to  include  a  second  Gaussian  term  in  the  model 
of  the  beam  profile,  but  there  is  a  parameter  that  must  be  tuned  to  adjust  the  relative 
weight  of  the  two  Gaussians.  This  was  discussed  at  PTCOG  with  Armin  Langenegger, 
one  of  Varian’s  physicists,  who  did  not  think  this  parameter  was  necessary  for  IBA  pencil 
beam.  This  second  Gaussian  was  added  for  the  Hitachi  beamline  at  MD  Anderson.  Next 
year,  Derek  will  continue  a  study  comparing  the  pencil  beam  profile  at  depth  between 
Eclipse  and  Monte  Carlo  calculations. 


Figure  BI-2.  Dose  distribution  for  a  prostate  case  calculated  with  Geant4 
Monte  Carlo  simulations.  Note  that  the  DVH  line  for  target  coverage  is 
much  worse  for  the  Monte  Carlo  calc  than  for  the  Eclipse  dose  calc.  This 
is  due  mostly  to  the  sample  variance  contribution  to  the  Monte  Carlo  dose 
uncertainty,  which  can  be  reduced  by  running  longer  simulations.  The 
median  dose  for  the  Monte  Carlo  calc  is  only  49  Gy  rather  than  the  5 1  Gy 
that  Eclipse  gives. 
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Figure  BI-3.  Dose  uncertainty  distribution  for  the  prostate  case  in  Figure 
1-2.  This  distribution  essentially  represents  the  error  bars  for  the  dose 
shown  in  Figure  1-2.  Note  that  the  uncertainty  shown  here  currently  only 
includes  the  contribution  of  random  sample  variance  inherent  from  the 
Monte  Carlo  method.  No  contribution  is  included  for  setup  uncertainties, 
motion,  target  delineation,  or  other  uncertainties.  The  white  areas  have 
uncertainty  less  than  4%.  The  top  right  panel  is  the  analogue  of  DVH:  an 
uncertainty  volume  histogram.  It  shows  that  the  target  dose  is  accurate  to 
between  2  and  5%.  The  relative  sample  variance  scales  as  l/sqrt(dose),  so 
the  most  uncertain  points  will  be  at  the  margin  of  the  PTV,  as  can  be  seen 
top  left  where  some  of  the  PTV  is  not  white,  i.e.  >  4%  sigma. 
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Figure  BI-4.  Distribution  of  the  dose  difference  between  Monte  Carlo  and  Eclipse.  A 
positive  difference  indicates  that  the  dose  from  Monte  Carlo  simulation  is  larger  than  the 
dose  from  Eclipse.  Red  areas  are  where  Monte  Carlo  is  5  Gy  (10  %  of  prescription) 
higher  than  Eclipse,  purple  areas  are  where  Monte  Carlo  is  5  Gy  lower  than  Eclipse. 
Simulation  is  showing  that  near  the  margin  of  the  PTV,  the  dose  is  generally  about  10% 
lower  than  Eclipse  (the  purple  border),  and  then  about  a  centimeter  beyond  the  margin  it 
is  10%  higher  (the  red/orange  halo).  We  believe  this  is  the  effect  of  a  dose  halo  that 
Monte  Carlo  is  getting  right,  but  Eclipse  is  failing  to  reproduce. _ 


We  can  also  do  the  forward  calc  for  double-scattering  and  uniform  scanning  plans.  Chris 
Ainsley  has  been  developing  the  simulation  code  for  the  universal  nozzle  since  he  began 
in  2006.  He  has  been  primarily  concerned  with  double-scattering  mode,  and  so  his 
simulation  code  includes  the  scatterers,  the  modulator  wheels,  and  more  recently  the 
MLC  and  patient-specific  compensator.  Derek  has  incorporated  the  voxelized  phantom 
based  on  patient  CT  data  into  the  Chris’s  universal  nozzle  code.  These  new  options  were 
added  to  the  macro  file  for  universal-nozzle  simulations: 


#===================================== 

#  detectorType: 

#  0  =  none 

#  1  =  readout  (default) 

#  2  =  voxelized  patient  (need  to  specify  CTDataPath) 

#  3  =  gas  chamber  layer 

#=========================== 


34 


/beamLine/Phantom/detectorType  2 

/beamLine/Phantom/CTDataPath  ct/datafiles/ 

/beamLine/Phantom/planFile 

ct/datafiles/PBS2/RN.  1.2.246. 352. 71. 5. 179413538 1. 322230.20100326105405. dcm 
/beamLine/Phantom/structureFile 

ct/datafiles/PBS2/RS.  1.2.246. 352.71. 4. 1794135381. 25548.20100326105458.dcm 
/beamLine/Phantom/bodyContour  BODY 

These  commands  setup  a  voxelized  geometry  based  on  the  dicom  CT  datafdes  in  the  path 
ct/datafiles/.  The  geometry  is  properly  oriented  in  the  beamline  based  on  the  isocenter, 
gantry  angle,  and  other  data  in  the  dicom  plan  file.  The  plan  file  also  contains  the  MLC 
shape  and  compensator  pattern.  The  bodyContour  is  just  used  to  set  the  material  to  air 
outside  the  specified  structure.  The  voxelized  CT  geometry  code  has  been  discussed  in 
previous  reports  and  will  not  be  repeated  here. 

Maura  Kirk  implemented  the  C++  class  to  write  out  dicom  RT  Dose  files.  This  means 
that  dose  scored  on  the  voxelized  patient  CT  geometry  can  be  written  to  a  dicom  file  that 
can  be  opened  by  a  dicom  viewer,  registered  with  the  CT  data,  with  dose  drawn  as 
colorwash,  manipulated  in  3D,  allowing  point  measurements,  DVHs  to  be  calculated,  and 
support  taking  differences  between  dose  maps  which  will  be  useful  as  we  start  to 
compare  Eclipse  calculated  dose  with  simulations.  We  are  using  Velocity  AI  to  open  the 
dicom  files  (CT  data  and  Monte-Carlo-calculated  RT  Dose  files).  Velocity  can  also  open 
dicom  structure  files,  and  plot  the  physician-drawn  contours  with  the  CT  and  dose  data. 
An  example  is  shown  in  Figure  1-1  where  a  Monte  Carlo  dose  calculation  is  being 
compared  to  the  Eclipse  calculation  in  Velocity.  The  plan  in  the  Figure  was  calculated 
with  Eclipse,  exported  as  a  DICOM  file  used  as  input  to  our  Monte  Carlo  code,  which 
outputs  a  DICOM  RT  Dose  file  that  can  be  read  and  displayed  with  Velocity. 

We  continue  to  refine  our  implementation  of  the  IBA  nozzles  in  simulation  to  reach 
better  agreement  with  data  collected  during  commissioning  of  the  proton  treatment 
rooms.  So  far  we  have  double-scattering  and  uniform  scattering  data  from  4  of  our 
treatment  rooms.  We  will  begin  accepting  the  last  room,  which  is  dedicated  PBS,  in  the 
summer  of  2011.  We  are  relying  on  IBA’s  specifications  for  spot-size  and  beam 
emittance  to  generate  realistic  beam  data  using  Monte  Carlo  that  can  be  used  to 
commission  Eclipse  for  research  purposes.  For  the  simulation  of  pencil  beam  delivery,  it 
was  learned  that  the  use  of  Fermi-Eyges  approximations  for  the  multiple  scattering 
contribution  is  not  adequate  to  determine  an  initial  phase  space  for  the  simulations.  We 
had  used  Fermi-Eyges  theory  to  quickly  compute  for  each  beam  energy  a  phase  space  for 
a  proton  beam  focused  to  a  specified  spot  size  at  isocenter  (see  Quarterly  Report  for  the 
1/1/09-3/31/09  period).  Fermi-Eyges  assumes  that  the  beam  profile  is  double-Gaussian 
and  remains  so  even  after  propagating  through  material.  The  multiple-scattering 
implementation  in  Geant4  uses  a  different  theory,  Lewis  theory,  that  is  exact  in  the  sense 
that  any  errors  in  the  moments  of  the  multiple-scattering  distributions  can  be  traced  back 
to  approximations  made  to  the  differential  cross  section.  Lewis  theory  computes  the 
lateral  moments  as  well  as  the  angular,  whereas  Fermi-Eyges  only  provides  the  angular 
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moments.  Lewis  theory  does  not  assume  that  the  beam  profile  is  Gaussian,  nor  does  a 
Gaussian  profile  necessarily  remain  Gaussian  after  multiple  scattering.  We  found  that  if 
we  use  Fermi-Eyges  to  estimate  the  multiple-scattering  to  determine  an  initial  phase 
space,  using  that  phase  space  for  Geant4  simulations  does  not  result  in  a  beam  of  the 
proper  spot  size  and  position  of  the  double  waist.  Derek  Dolney  wrote  some  simulation 
code  to  track  an  infinitesimal,  monoenergetic  pencil  beam  through  air  and  fit  a  function 
to  the  multiple-scattering  contribution  as  a  function  of  depth  in  air.  This  fit  is  now  used 
instead  of  Fermi-Eyges  to  get  the  multiple-scattering  right  for  Geant4  simulations  in 
order  to  produce  a  beam  of  the  right  transverse  dimensions,  within  0.4%  (see  Figures  1-5 


o  0.5  1  1.5  2 

depth  in  air,  m 


Figure  BI-5.  Contribution  to  the  beam  sigma  due  to  multiple  scattering  as  a 
function  of  depth  in  air.  The  lines  represent  the  Fermi-Eyges  calculation,  the 
spots  are  Geant4  simulations  which  use  Lewis  scattering  theory.  Fermi- 
Eyges  overestimates  the  scattering  in  air,  so  that  if  a  phase  space  is 
constructed  using  the  Fermi-Eyges  estimate,  the  double  waist  of  the 
simulated  beam  will  be  deeper  than  was  desired,  as  can  be  seen  in  Figure  I- 
6. 
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Figure  BI-6.  Spot  size  near  isocenter  from  simulations  using  a  primary 
phase  space  calculated  with  Fermi-Eyges  (lower  line)  vs.  using  a  fit  to  the 
multiple-scattering-in-air  contribution  from  Geant4  simulations.  The  target 
spot  size  at  isocenter  is  indicated  by  the  X.  Fermi-Eyges  overestimates  the 
multiple  scattering  in  air,  so  that  when  a  phase  space  is  calculated  with 
Fermi-Eyges,  the  simulated  waist  lies  at  a  greater  depth  in  air. 

Derek  has  implemented  the  geometry  for  the  IBA  dedicated  nozzle  in  Geant4  and 
determined  a  phase  space  for  primaries.  This  nozzle  results  in  a  smaller  spot  size  because 
the  entire  nozzle  is  evacuated  (Figure  1-7).  Initial  phase  space  was  found,  using  the 
procedure  described  above,  to  fit  the  IBA  specification  spot  size  at  isocenter. 
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Figure  BI-7.  Spot  size  for  the  two  IBA  nozzles  as  a  function  of  range  in 
water.  Also  shown  is  the  spot  size  when  the  7.5  cm  acrylic  range  shifter  for 
treating  shallow  depths  is  inserted  into  the  beam  line. 


Derek  implemented  in  simulation  the  7.5  cm  range  shifter  solution  that  IBA  is  using  to 
treat  shallow  targets  with  pencils  beams.  Derek  generated  commissioning  data  and 
learned  how  to  commission  Eclipse  with  range  shifter  data  and  to  setup  and  calculate 
treatment  plans  that  use  the  range  shifter.  Our  treatment  planning  system  is 
commissioned  with  Monte  Carlo  data  for  the  universal  nozzle  and  dedicated  nozzle, 
including  the  range  shifter,  and  we  have  begun  to  perform  studies  that  compare  PBS 
treatment  plans  with  double-scattering  plans  and  with  IMRT  plans.  Several  abstracts  were 
submitted  to  PTCOG  and  ASTRO  for  pancreatic,  head  and  neck,  breast,  and  rectal  plan 
comparison  studies. 

Next  year  will  be  a  very  exciting  time,  as  we  will  finally  have  PBS  data  to  validate  the 
Monte  Carlo  simulations.  We  expect  some  refinements  to  the  nozzle  will  be  necessary. 
Certainly  it  will  be  necessary  to  generate  new  initial  phase  space  data  since  IBA  is  saying 
that  beam  is  not  perfectly  focused  at  isocenter  but  as  much  as  20  cm  upstream.  This  can 
be  accommodated  easily  by  the  procedure  described  above  and  in  previous  reports  once 
the  data  is  available.  Perhaps  we  will  also  find  that  some  of  the  simplifications  we  have 
made  (geometric  rather  than  true  magnetic  beam  steering,  for  example)  are  too  simplistic 
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and  need  to  be  improved.  Until  the  data  comes,  we  still  have  differences  between  Eclipse 
dose  calc  and  Monte  Carlo  that  need  to  be  understood. 

II.Development  of  a  Micromegas  layer  for  proton  dosimetry  in  2D 

Derek  Dolney  has  worked  on  the  design  and  simulation  of  a  system  of  detectors  for  three- 
dimensional  dosimetry  in  proton  beams.  The  proposed  detector  system  is  based  on  a  type 
of  ionization  chamber  called  Micromegas.  Micromegas  uses  a  micromesh  screen  held  at 
high  voltage  to  achieve  signal  amplification  in  a  relatively  thin  gas  chamber  (see  Figures 
II-2,  II-3,  and  II-4).  The  ionization  current  can  be  collected  on  a  pixelated  anode  for  2D 
dose  resolution.  Three-dimensional  resolution  is  obtained  by  stacking  multiple  2D  planes 
back-to-back.  Very  large  systems  of  Micromegas  detectors  are  used  in  high-energy 
physics  experiments  for  high-resolution  particle  tracking,  but  the  dose  rate  is  smaller  than 
the  dose  rate  in  the  proton  Bragg  .  It  should  be  possible  to  achieve  a  Micromegas  design 
for  a  proton  therapy  beam  that  achieves  the  same  high-resolution,  but  does  not  saturate  in 
the  Bragg  peak,  by  reducing  the  gain.  The  signal  gain  is  a  function  of  operating 
potentials,  gap  thicknesses,  and  gas  mixture  and  pressure.  There  is  already  some  data  in 
the  literature  about  the  dependence  of  the  gain  on  these  parameters,  and  it  is  expected, 
based  on  the  energy  deposit  in  the  detector  obtained  from  Monte  Carlo  simulations  (see 
Figure  II- 10),  that  a  successful  design  can  be  achieved.  Eventually,  multiple  Micromegas 
layers  stacked  back-to-back  would  provide  3D  dose  resolution. 
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Figure  BII-1.  Schematic  of  a  Micromegas  detector  for  2D  proton  dosimetry  (left). 
Ionization  is  produced  in  the  drift  gap  by  the  incoming  proton  beam.  Those  ions  drift 
towards  the  micromesh.  The  large  electric  field  in  the  vicinity  of  the  micromeshand  a 
system  for  3D  dosimetry  built  with  2D  Micromegas  layers  interleaved  with  layers  of 
solid  water  material  (right). 


A  prototype  Micromegas  chamber  has  been  assembled  in  the  Department  of  Physics  at  the 
University  of  Pennsylvania  (Figures  II-2,  II-3,  II-4).  Currently  the  anode  is  a  single  channel. 
This  prototype  is  current  in  a  testing  phase  at  in  the  Physics  Department. 
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Figure  BII-2.  Pictures  of  the  Micromegas  detector  prototype  assembled  in  the 
Department  of  Physics  at  the  University  of  Pennsylvania.  Left:  the  gray  rings  seal  the 
electrodes  and  micromesh  with  O-rings.  Gap  spacings  can  be  changed  by  mounting  the 
electrodes  and/or  the  mesh  in  different  slots  in  the  gray  mounting  rings.  A  wire  for  one 
of  the  high-voltage  electrodes  is  visible.  Right:  with  one  of  the  electrodes  removed,  the 
micromesh  is  visible  inside  the  chamber. 


Figure  BII-3.  Pictures  of  the  Micromegas  detector  prototype  assembled  in  the 
Department  of  Physics  at  the  University  of  Pennsylvania.  Left:  the  inner  blue  ring 
maintains  tension  on  the  micromesh  to  prevent  it  from  deforming  in  the  electric  field. 
Right:  another  view  of  the  micromesh.  Also  visible  are  the  additional  slots  in  the  gray 
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spacing  rings  that  allow  the  gap  thicknesses  to  be  varied. 


Figure  BII-4.  An  earlier,  unsatisfactory 
attempt  to  prevent  deformation  of  the 
micromesh  in  the  electric  field.  The  outer 
blue  ring  was  tightened  to  the  gray  spacing 
rings  to  secure  the  mesh.  This  caused  the 
mesh  to  buckle,  rather  than  pull  flat,  and 
was  abandoned  in  favor  of  the  ring  in 
Figure  II-3  (left).  One  of  the  sealing  O- 
rings  (black)  can  be  seen  in  this  photo. 


Two  experiments  are  underway  to  calibrate  the  gas  gain  of  the  detector  system:  1.  An  Fe- 
55  source  is  being  used  to  calibrate  the  gas  gain,  and  2.  High-energy  electron  beams 
produced  by  the  linear  accelerators  in  the  Department  of  Radiation  Oncology  are 
minimum  ionizing,  or  nearly  so,  and  can  be  used  to  obtain  a  gas  gain  curve.  We  need  to 
characterize  the  gas  gain  for  Micromegas  in  a  proton  beam  so  that  design  choices  such  as 
operating  potentials,  gas  composition  and  pressure,  and  gap  thicknesses  can  be  made.  We 
can  use  simulations  and  theory  to  relate  the  detector  response  to  electrons  with  the 
response  to  protons,  though  we  also  verify  with  measurements.  A  low-voltage  gain  curve 
is  shown  in  Figure  11-5.  We  used  an  electrometer  to  provide  the  bias  voltage.  The 
electrometer  only  goes  up  to  500  V.  We  need  higher  voltage  to  find  the  proportional 
region  of  the  response  curve.  We  have  external  power  supplies  that  will  go  to  5  kV,  and 
obtained  a  special  adapter  from  IBA  Dosimetry  in  order  to  operate  our  electrometer  in  a 
grounded  input  mode  that  will  work  with  an  externally  supplied  bias.  Derek  found  the 
response  curve  to  be  essentially  flat  to  4500  V.  Above  that  level  the  prototype  began 
arcing  and  could  not  be  used  reliably.  A  second  prototype  was  built  (see  below)  that  has 
some  engineering  improvements  that  should  allow  stable  operation  at  higher  voltages. 
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Bias  Voltage 

Figure  BII-5.  Gain  curve  for  the  prototype  Micromegas  chamber  in  a  20 
MeV  electron  beam. 


Some  data  has  been  collected  with  the  prototype  Micromegas  detector  in  the  proton 
beam,  and  in  electron  beams.  Initially,  a  brief  look  at  the  reproducibility  of  the  prototype 
in  the  proton  beam  was  considered.  Those  results  are  given  in  Figure  11-13.  The  response 
of  the  detector  slowly  decreases  over  time.  We  believe  this  is  due  to  a  space-charge 
effect:  charge  is  collecting  on  non-conducting  surfaces  in  the  prototype  assembly,  such  as 
the  o-rings  and  plastic  chamber  housings,  thus  lowering  the  local  electric  field.  The 
problem  appears  to  be  reduced  by  operating  at  a  higher  bias  (300V).  The  next  prototype 
will  have  guard  ring  to  better  shape  the  field  around  the  edge  of  the  detector  so  that 
charge  is  collected  on  the  electrodes  rather  than  accumulating  on  insulators. 
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Figure  BII-6.  Measurements  with  the  prototype  detector  at  two  bias  voltages,  attempting 
to  establish  measurement  reproducibility  in  terms  of  charge  collected  per  proton  beam 
monitor  unit  (MU).  At  150  V,  the  detector  response  decreases  with  time,  suggesting  that 
we  are  collecting  charge  on  the  non-conductive  materials  around  the  edges  of  the  detector 
active  volume,  the  so-called  “space-charge  effect”.  Operating  at  300  V  seems  to  give 
better  reproducibility.  Our  second  prototype  will  have  a  guard  ring  to  more  efficiently 
collect  charge  near  the  edges. _ 


A  proton  depth  dose  curve  was  measured  with  the  prototype  Micromegas  detector  (Figure 
II-7)  and  compared  with  the  depth  dose  measured  with  the  commercial  PTW  "Bragg  Peak 
Chamber".  The  two  curves  were  normalized  in  the  entrance  region  (at  about  8  cm)  and 
shifted  relative  to  each  other  in  the  depth  direction,  because  the  range  is  not  well 
determined  for  this  type  of  delivery.  The  aim  was  to  verify  that  a  Micromegas  design 
does  not  saturate  in  the  high  dose  rate  within  the  Bragg  peak.  The  data  plotted  in  Figure 
II-7  suggests  that  we  are  on  the  right  track,  but  a  more  careful  measurement  is  needed. 

We  need  to  take  more  data  points  in  the  vicinity  of  the  Bragg  peak  with  the  prototype. 
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Figure  BII-7.  A  measurement  of  a  proton  prisitine  Bragg  peak  with  the 
prototype  chamber.  A  measurement  taken  with  the  PTW  Bragg  Peak 
chamber  is  shown  for  comparison.  The  curves  are  normalized  in  the 
entrance  region  at  8  cm.  The  two  measurements  were  shifted  in  the  depth 
direction  due  to  the  range  uncertainties  associated  with  this  type  of  beam 
delivery. 


A  second-generation  prototype  has  been  assembled  by  Professor  Robert  Hollebeek's 
group  in  the  Physics  Department.  This  prototype  has  several  improvements  over  the  first 
one.  The  most  significant  difference  is  that  the  cathode  of  the  new  prototype  is  segmented 
into  16  strips  so  that  we  can  begin  to  make  measurements  to  determine  spatial  resolution 
for  this  Micromegas  application.  This  second  pilot  device  also  has  several  high  voltage 
and  mechanical  improvements.  This  second  device  will  also  allow  us  to  have  one  device 
we  can  put  in  the  beam  and  still  keep  one  on  the  test  stand.  Eventually  we  will  make  a 
prototype  with  pads.  The  electronic  readout  for  the  strips  is  under  construction  and  there 
is  an  un-segmented  plane  that  can  be  used  in  the  meantime.  The  whole  assembly  can  be 
taken  apart  in  a  few  minutes  and  the  ring  arrangement  allows  us  to  quickly  change  the 
spacing  of  the  cathodes  and  screen. 

The  sequence  of  pictures  below  shows  the  layers  in  this  prototype  proceeding  from  top  to 
bottom.  Each  picture  has  one  additional  layer  removed. 
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1.  The  top  bolt  flange  keeps  pressure  on  the  top  ground  and  cathode  layers.  This  plate  is 
1/32”  double  sided  G10.  The  top  side  is  a  ground  shield  and  the  lower  side  is  the  top 
cathode. 
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3.  Gas  inlet  ring  and  o-ring  seal.  This  version  has  a  gas  inlet  ring  which  eliminates 
positive  charge  buildup  on  the  gas  tubes  caused  by  having  the  gas  feeds  enter  at 
the  top  of  the  chamber.  Eventually  the  charge  buildup  leads  to  a  corona  discharge 
to  the  cathode  or  ground.  The  gas  inlet  holes  are  on  the  sides  and  not  visible  in 
this  view. 
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4.  Screen  ring.  The  screen  now  has  a  machined  guard  ring  to  prevent  corona  at  the 
edges  of  the  screen  where  there  are  14  micron  wire  tips  that  have  to  be  terminated 
in  a  smooth  surface.  The  stainless  screen  is  soldered  to  the  copper  guard  ring. 

The  removable  screen  assembly  sits  in  a  groove  in  the  ring. 
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5.  Separator  layer  used  to  adjust  the  screen  to  cathode  distance  and  the  segmented 
cathode  layer.  The  cathode  layer  has  ten  segments  on  the  right  (2  guard  regions 
and  8  strips)  and  ten  segments  on  the  left.  The  left  right  arrangement  will  allow  us 
to  measure  beam  centering  and  jitter  and  the  strips  will  allow  us  to  test  the 
position  resolution  of  the  cathode.  The  strip  pitch  is  3  strips  per  centimeter  so  we 
estimate  the  resolution  to  be  better  than  a  tenth  of  a  millimeter.  The  serial  readout 
will  have  a  time  resolution  of  about  200kHz  per  segment  giving  a  complete 
readout  every  100  microseconds. 
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6.  Side  view  (unassembled).  The  layers  are: 
Top  flange 
Top  cathode 
Gas  feed  ring 
Screen  layer 
Spacer 

Segmented  cathode 
Bottom  flange 


The  Poisson  solver  software  Poisson  Superfish  was  used  to  obtain  plots  of  the  field 
strength,  given  the  detector  geometry  and  potentials.  Figure  II-8  shows  the  high  field  in 
the  vicinity  of  the  micromesh,  where  signal  amplification  will  occur.  Figure  11-16  plots 
the  field  strength  along  the  central  axis  of  the  detector.  The  signal  gain  in  the  detector  is 
related  to  the  field  strength. 
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Figure  BII-8.  Electromagnetic  field  line  plot  obtained  with  the  Poisson  Superfish 
software  using  the  prototype  detector  geometry.  The  mesh  is  12  um  wires  at  0.025  cm 
above  the  anode.  The  signal  amplification  is  highest  near  the  micromesh,  where  the  field 
lines  are  most  dense.  Additional  amplification  is  achieved  in  the  amplification  gap. 


Electric  field  data  from  file  PLATE . AH 
Psoblem  tiae  1:  Me:>h  wise  fields 


Y  ( cm ) 


Figure  BII-9.  Field  strength  along  the  central  axis  from  Poisson  Superfish  (see  Figure 
11-15)  showing  the  high  field  in  the  vicinity  of  the  micromesh  at  0.025  cm. 
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Figure  BII-10.  Energy  deposit  in  the  detector  for  different  gas  mixtures  from  Geant4 
Monte  Carlo  simulation.  The  gas  pressure  is  1  atm.  The  beam  current  is  assumed  to  be  3 
nA,  which  is  typical  for  double-scattering  and  scanned  pencil  beam  deliveries. 


We  will  develop  this  Micromegas  design  very  actively  over  the  next  year.  We  are  working 
with  CERN  to  obtain  two  small  Micromegas  prototypes.  These  prototypes  will  have  narrow 
amplification  gaps,  allowing  higher  field  strength  and  gas  gain.  The  CERN  group  are  the 
experts  in  manufacturing  Micromegas  layers  in  bulk.  We  will  measure  depth  dose  curves  and 
beam  profiles  with  the  2D  layers.  Our  goal  by  the  end  of  May  2012  is  to  have  a  system  of  4 
Micromegas  layers  of  at  least  20  cm  x  20  cm  transverse  size  that  can  measure  a  20  cm  x  20 
cm  x  4  cm  segment  of  any  dose  distribution  (generated  with  Eclipse  or  other  planning 
system). 
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C.  Image-Guided  and  Adaptive  Proton  Radiotherapy 


At  this  time,  we  have  requested  a  no-cost  extension  with  a  new  Scope-of-Work  to  perform  this 
adaptive  study.  We  are  still  waiting  for  a  response  to  this  request. 


Key  Research  Accomplishments 

The  major  accomplishment  of  the  past  year  has  been  the  successful  testing  of  the  multileaf 
collimator.  While  there  is  still  one  component  to  be  delivered  and  more  testing  to  be  performed  it 
is  clear  that  the  MLC  will  be  a  more  efficient  method  to  collimate  the  proton  beam  for  the 
majority  of  fields  than  using  a  machined  aperture. 


A  second  accomplishment  has  been  the  robust  growth  in  clinical  treatment  protocols  for  our 
patients.  There  are  currently  15  IRB-approved  proton  therapy  protocols  accruing  patients  at  our 
center.  A  total  of  189  adults  and  44  CHOP  patients  have  completed  proton  treatment  as 
planned.  The  total  number  of  completed  individual  proton  treatments  is  233.  At  present,  on 
treatment  there  are  48  adults  and  9  CHOP  patients  for  a  total  of  57  patients. 

Reportable  Outcomes 

The  following  abstracts  based  on  work  performed  on  this  project  have  been  accepted 
during  the  past  year  at  scientific  meetings: 

1.  Optimizing  Monte  Carlo  Parameters  to  fit  Proton-beam  Measurements  with  a  Tungsten 
MLC,  PTCOG  2011 

2.  Monte  Carlo  Studies  of  Double  Scattering  Proton  Beam  Characteristics  with  an  MLC, 
DYAAPM  Young  Investigator  Symposium 


Conclusions 

This  report  documents  the  work  that  has  been  accomplished  during  the  sixth  year  of  the 
project  to  design  an  MLC  for  proton  radiotherapy,  the  fifth  year  of  work  on  the  scanned 
beam  project,  and  the  fourth  year  of  work  on  the  image-guided  proton  therapy  project.  It 
concentrates  on  the  past  quarter  since  reports  on  the  other  quarters  already  have  detailed 
those  efforts. 
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